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Abstract 

Using a new type of simulation algorithm for the standard SU(3) lattice gauge theory that 
yields results with unprecedented precision, we confirm the presence of a j/r correction 
to the static quark potential at large distances r, with a coefficient 7 as predicted by the 
bosonic string theory. In both three and four dimensions, the transition from perturbative 
to string behaviour is evident from the data and takes place at surprisingly small distances. 



1. Introduction 

The idea that the pure SU(3) gauge theory (and all other non-abelian gauge theories) 
might be closely related to some kind of string theory has been around for a long 
time. A comparatively direct connection derives from the observation that the gauge 
field generated by two widely separated static quarks appears to be squeezed into a 
flux tube. Such flux tubes may be expected to behave like strings with fixed ends, 
at least in the limit where they are very much longer than wide. 
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On a more formal level it has been suggested that the expectation values of large 
Wilson loops are matched by the corresponding amplitudes of an effective bosonic 
string theory [1]. Assuming this to be the case, the static quark potential V(r) can 
be shown to have an asymptotic expansion 

V(r) = ar + l i + 1 /r + 0(l/r 2 ), 7 = -^(d-2), (1.1) 

at large distances r, where a denotes the string tension, jjl a regularization-dependent 
mass and d the dimension of space-time [2,3]. The 7/r correction in this formula is a 
quantum effect that is characteristic of the relativistic bosonic string. It is universal 
in the sense that the value of 7 is the same for a large class of string actions. 

In principle the validity of eq. (1.1) can be checked by putting the theory on a 
lattice and by calculating V(r) through numerical simulation. The statistical and 
systematic errors in these computations are, however, rapidly increasing at large 
distances, and it is then practically impossible to clearly separate the 7/r correction 
from the other terms (if the standard simulation algorithms are used). As a result, 
we would find it difficult to claim, on the basis of the published simulation data 
alone, that the expansion holds with the quoted value of 7. 

There is nevertheless some support from lattice gauge theory for the string model. 
It has been noted, for example, that the data for the potential at distances between 
0.4 and 0.8 fm agree with eq. (1.1) within small errors [4]. Another numerical study 
that we wish to mention here is the work of Caselle et al. [5] on the confinement phase 
of the Z2 gauge theory in three dimensions. Highly efficient simulation techniques are 
available in this case, and using these it was possible to confirm that the expectation 
values of large rectangular Wilson loops are matched by the string theory amplitudes 
to a precision where the subleading string effects can be clearly identified. 

Last year we proposed a new type of simulation algorithm for the SU(3) gauge 
theory that leads to an exponential reduction of the statistical errors in calculations 
of the Polyakov loop correlation function [6]. Similar to the expectation values of 
Wilson loops, the latter can serve as a probe for string effects [7,8], and the onset 
of string behaviour may in fact be easier to observe in this way [9]. In the present 
paper, though, the correlation function only appears at an intermediate stage. The 
computations are well adapted to the new algorithm, and the results that we obtain 
for the force V'(r) and the second derivative V"(r) are sufficiently precise for the 
question raised above to be addressed. 
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2. Polyakov loop correlation function 

In this section we recall some well-known facts about the Polyakov loop correlation 
function in the SU(3) gauge theory. The notation is the same as in ref. [6], except 
that we now consider the theory also in three dimensions. 

2.1 Definition 

We set up the lattice theory in the standard manner on a 4-dimensional hypercubic 
or 3-dimensional cubic lattice with spacing a, time-like extent T and spatial size L. 
Periodic boundary conditions are imposed in all directions, and for the gauge field 
action we take the Wilson plaquette action S[U] with bare gauge coupling g [10]. 
For any given gauge field configuration U(x,/jl), the Polyakov loop 

P(x) = tr {U(x, fi)U(x + a£, fi) . . . U(x + (T - a)£, »)}^ =0 (2.1) 

is the trace of the Wilson line that passes through x and that winds once around the 
world along the time axis (/t denotes the unit vector in direction fi). We are then 
interested in the correlation function 

{P(x)*P(y)) = 1 f D[U]P(xyP(y)e- s M, B[U] =Y[dU(x,»), (2.2) 

where Z is a normalization factor such that (1) = 1 and dU(x,fi) stands for the 
normalized invariant measure on SU(3). 

2.2 Spectral representation 

An important property of the Polyakov loop correlation function is that it can be 
written as a ratio of partition functions. The transfer matrix formalism must be 
invoked to show this, but the derivation is then fairly straightforward (for an intro- 
duction to the subject, see refs. [11-14] for example). 

First recall that the quantum mechanical states of the theory are represented by 
wave functions on the space of all spatial lattice gauge fields U at time xo = 0. The 
transfer matrix in the temporal gauge, T = e _Ha , is a bounded hermitian operator 
that acts on these functions. In this framework the gauge-invariant wave functions 
describe the physical states in the vacuum sector, and the wave functions with two 
colour indices that transform according to 

MU A ] = A(x) a7 [A(yM*Yv[£/] (2.3) 
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under gauge transformations A (where U denotes the gauge transform of U) de- 
scribe states with a static quark at the point x and a static antiquark at y. 

If we introduce the projector P to the subspace of gauge-invariant wave functions, 
and similarly the projector P303* (x, y) to the subspace of wave functions (2.3), the 
partition functions in these sectors are given by 



Z = Tr{Pe- HT }, 

Zsm* (x, y) = \ Tr {P 3 ®3* (x, y) e" HT } . 



(2.4) 
(2.5) 



A conventional normalization factor has been included here to compensate for the 
trivial degeneracy of the states in the multiplets (2.3). The partition function (2.5) 
is thus effectively a sum over multiplets. 

By inserting the explicit expression for the transfer matrix, both partition func- 
tions can be rewritten in the form of functional integrals on the d-dimensional lattice. 
Using an SU(3) character relation, it is then possible to establish the identity 

(P(x)*P(y)) = z ^f' y) . (2.6) 

The spectrum of the transfer matrix is purely discrete in finite volume, and since 
the ground state in the vacuum sector is non-degenerate [12,13], an immediate con- 
sequence of this representation is that 



(P(x)*P(y)) = J2^ne- 



E n T 



71 = 



(2.7) 



for some positive energies E n and integral weights w n . 
2.3 Relation to the potential V(r) 

We shall only consider the potential along the lattice axes and thus set y = x + r/t, 
fj, = 1, in the following. If we arrange the energy values in eq. (2.7) in ascending 
order, it is clear that E is the lowest energy (above the ground state energy) in the 
sector with the static quarks. In other words, E coincides with the static quark 
potential V(r). 

The associated eigenstate can be shown to be non-degenerate at strong coupling, 
and in our numerical studies of the Polyakov loop correlation function we have always 
found that wq = 1. It is usually fairly easy to obtain this result, since wq is known 
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to be a positive integer. We thus conclude that 



V{r) = ~\n{P{xyP{y)) + e 



(2.8) 



where 




1 



-AET 



+ ...} 



(2.9) 



In particular, no fits are required to extract the potential from the Polyakov loop 
correlation function once the latter has been computed at large values of T, where 
the error e is negligible. 



3. Bosonic string theory 

The following paragraphs contain a brief description of the string theory represen- 
tation of the Polyakov loop correlation function. It is included here mainly for the 
reader's convenience, but also to prepare the ground for the analysis of the simulation 
data. 

3.1 Effective action 

In the bosonic string theory the Polyakov loop correlation function is formally given 
by a functional integral over all two-dimensional surfaces (world-sheets) that are 
bounded by the loops [1]. When both T and r are large with respect to the scale set 
by the string tension a, the integral is dominated by the surface with the minimal 
area and can be calculated by expanding the integrand around this surface [2]. 

The important degree of freedom in this situation is the deviation of the fluctuating 
surface from the minimal surface in the directions orthogonal to the latter (fig. 1). If 
we parameterize the minimal surface through two cartesian coordinates z a (a = 0, 1) 
in the obvious way, the deviation is a vector field 



with d— 2 non-zero transverse components that satisfy Dirichlet boundary conditions 
at z\ = and z\ = r. The associated effective action may depend on how precisely 



h{z) = (0,0,/i2(z),...,/id-i(z)) 



(3.1) 
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Fig. 1. In the effective string theory, the Polyakov loop correlation function is ob- 
tained by expanding the associated string partition function about the minimal surface 
(shaded area) spanned by the loops. 

the underlying bosonic string theory was set up, but for symmetry reasons it has to 
be of the form 



where the ellipses stand for higher-dimensional interaction terms that are Lorentz- 
invariant polynomials in the derivatives of h [3]. In this context the field h is counted 
as dimensionless and it is, therefore, natural to assume from the beginning that the 
displacement of the fluctuating surface in physical units is h/\fa rather than h. 

From a purely technical point of view, the effective string theory outlined above 
is a two-dimensional quantum field theory on the world-sheet parameter space with 
specified boundary conditions. Since the coefficients of the possible interaction terms 
are dimensionful, the perturbation expansion generated by the effective action is non- 
renormalizable, and counterterms of increasing dimension thus need to be included 
to cancel the divergent parts of the loop diagrams. 

3.2 Leading- order approximation 

If we neglect the interaction terms in the effective action (3.2), the functional integral 
over all vector fields h(z) can be evaluated exactly. This leads to the expression 



where A denotes the laplacian on a two-dimensional cylinder of height r and cir- 
cumference T with Dirichlet boundary conditions. The mass jj, that appears in this 
formula is a free parameter of the string model whose value has no particular physi- 
cal meaning. In the gauge theory there is in fact a similar normalization ambiguity, 
once the need for renormalization of the Polyakov loop correlation function is taken 
into account [15,16]. 




(3.2) 



(P(x)*P(y)) = e 



'-" T [det(-A)]-*< d - 2 > , 



(3.3) 
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The determinant of the laplacian can be computed explicitly using the £-function 
finite-part prescription or one of the standard regularization methods. Such calcula- 
tions have been performed in refs. [17,7,5], for example, and also in many papers on 
conformal field theory and fundamental string theory (see refs. [18,19] for a review 
and a list of references) . In terms of the Dedekind r\- function, 

oo 
n=l 

the result is 

det(-A) = ^) 2 , q = e~^ r , (3.5) 

from which we immediately deduce that the string amplitude on the right-hand side 
of eq. (3.3) admits an expansion of the form (2.7), with energies and weights 

E n = ar + fi + ^{-±(d-2)+n}, (3.6) 



w 



1, Wl = d-2, w 2 = |(d-2)(d + l), etc. (3.7) 



In particular, to this order of the effective string theory, the static quark potential 
is given by eq. (1.1) (without any corrections proportional to higher powers of 1/r) 
and the energy gap AE is equal to n/r. 

3.3 Higher orders & renormalization 

Quantum field theories on space-time manifolds with boundaries require, in general, 
the inclusion of terms in the action that are localized at the boundary [20,21]. In 
the present case the simplest term of this type is 

S x = \b\ dz {{dihdrh)^ + (drhdih)^} , (3.8) 
J o 

where 6 is a parameter of dimension [length] . The effect of Si on the Polyakov loop 
correlation function can be easily worked out to first order in b using the £-function 
finite-part prescription, for example. For the static potential and the energy gap the 
calculation yields 

V(r) = or + /x - ^ (d - 2) (1 + b/r) , (3.9) 
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AE = — (1 + b/r) . 



(3.10) 



As expected these corrections are negligible compared to the leading-order terms at 
sufficiently large values of r. 

It can be shown that the boundary action S\ is the only term with a coefficient of 
dimension 1 that can occur in the effective action. There are just a few more terms 
with coefficients of dimension [length] 2 , such as 

S 2 = \c2 ! [ dz d Zl (d a hd a h) (d b hd b h) . (3.11) 
Jo Jo 

They give rise to corrections to the static quark potential proportional to 1/r 3 , i.e. to 
terms that are even more suppressed at large distances. 

A more careful discussion of the higher-order contributions must take into account 
the fact that the theory needs to be regularized. This introduces another scale, the 
ultra-violet cutoff A, and corrections proportional to (c2/V 3 )(Ar) p , for example, must 
then be expected to arise. However, as in any other non-renormalizable field theory, 
it must be possible to cancel the divergent contributions by a renormalization of a, 
\i and the coefficients of the interaction terms. 

For illustration we note that a term proportional to c 2 A/r 2 actually does occur 
if a lattice regularization is used. In this case the inclusion of the boundary term 
Si is thus required to be able to renormalize the theory. The bottom line is then 
that equations like (3.9) and (3.10) hold to all orders of the renormalized effective 
theory, up to corrections proportional to higher powers of 1/r. 



4. Numerical evaluation of V(r) and its derivatives 



Using the simulation algorithm of ref . [6] , the Polyakov loop correlation function can 
be computed at values of T and r that previously appeared to be inaccessible. In 
particular, when T is increased at fixed r, the computer time required for a specified 
relative statistical error is only growing approximately like T 3 , while the signal is 
exponentially decreasing. This enables us to obtain the static quark potential V(r) 
via eq. (2.8) with small errors. 

The effective string theory predicts that the force F(r) = V'{r) decreases towards 
the string tension a at large distances r with a rate proportional to 1/r 2 . In this 
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paper our principal goal is to find out whether the numerical data are consistent 
with this and whether the slope 



c(r) = \r A F'{r) (4.1) 

converges to the value of the universal coefficient 7 in eq. (1.1). We are, therefore, 
more interested in the derivatives of the static potential than in the potential itself. 

4-1 Lattice definition of F(r) and c(r) 

There is no unique way to "differentiate" lattice functions, but a sensible prescription 
should evidently be such that the correct continuum limit is obtained and that no 
artificially large lattice effects are introduced. The definition of F(r) proposed by 
Sommer [22] fulfils this condition and is easily extended to c(r). We are thus led to 
define 

F(f) = {V{r)-V{r-a)}/a, (4.2) 
c(f) = \? 3 {V(r + a) + V(r -a)- 2V{r)}/a 2 , (4.3) 

where f(r) and r(r) satisfy 

f = r -\a + 0(a 2 ), (4.4) 
r = r + 0(a 2 ). (4.5) 

The precise choice of these functions takes into account the approximate shape of 
the potential, so that an enhancement of lattice effects is avoided, particularly at 
small distances [22,4]. Explicitly we set 

r 2 ~ d = 2n(d - 2){G(r -a)- G(r)}/a, (4.6) 
f l ~ d = 2n{G{r + a) + G(r - a) - 2G(r)}/a 2 , (4.7) 

where G(r) denotes the value of the Green function of the lattice laplacian in d — 1 
dimensions at the point (r, 0, . . . , 0) (for a detailed discussion of the Green function 
in two and three dimensions, see refs. [23,4] for example). 
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4-2 Systematic errors 

Lattice spacing effects, the omission of the contribution of the higher-energy states 
in eq. (2.8) and finite-volume effects are sources of systematic uncertainties in our 
calculations. As usual the lattice effects can be controlled by performing simulations 
at several values of the lattice spacing, and we shall return to this issue in sect. 5 
when we discuss the simulation results. 

To estimate the error e in eq. (2.8), some information on the energy values E n 
and the associated weights w n is needed. Accurate data for the lowest excited levels 
were first obtained by Michael and Perantonis [24], and more precise and extensive 
studies have later been reported by Juge, Kuti and Morningstar [25,26]. The results 
show that, in four dimensions, the gap AE is about 2.2/r at the Sommer reference 
scale r = r . It then approaches ir/r around r = 3r and stays above this value at 
still larger distances. The associated weight w\ is 2 as in the effective string theory. 
Using this and the fact that the higher levels are well separated from Ei, the error 
in eq. (2.8) can be estimated straightforwardly. 

There are unfortunately no results on the spectrum of the excited levels in three 
dimensions, and we can only refer to the analogous case of the SU(2) theory, where 
the situation does not appear to be very different from the one described above [26]. 
We have, therefore, estimated the error in eq. (2.8) in the same way, with w\ = 1, 
and checked the estimate by performing simulations at various values of T. 

Finite-volume effects can arise when the colour flux tube is squeezed in the trans- 
verse directions by the finite extent of the lattice. In all our calculations the spatial 
lattice size L was at least 20 times larger than the bulk correlation length [27,28], 
and these effects were thus suppressed by many orders of magnitude. Another sys- 
tematic effect can be traced back to the presence of states where the flux tube is 
winding "around the world". The energy of these states is at least V(L — r), and 
their contribution to the error in eq. (2.8) will therefore be small if L and r are such 
that V(L — r) — V(r) is large. 

4-3 Statistical errors 

The parameters of the simulation algorithm of ref. [6] can be chosen so that subse- 
quent "measurements" of the Polyakov loop correlation function may be assumed to 
be statistically independent. We have done so and later confirmed the independence 
of the data by determining the integrated autocorrelation time. The standard single- 
elimination jackknife analysis was then applied to calculate the statistical errors. 

An important property of the data series generated in this way is that the values 
of the static potential at different distances turn out to be very strongly correlated. 
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Table 1. Simulation parameters 



d 


P 


r /a 


a [fm] 




T/a 


L/a 


4 


5.70 


2.93 


0.171 


7 


24 


18 


4 


5.85 


4.09 


0.122 


9 


36 


24 


4 


6.00 


5.39 


0.093 


12 


48 


30 


3 


11.0 


3.30 


0.152 


9 


32 


24 


3 


15.0 


4.83 


0.104 


12 


48 


32 


3 


20.0 


6.71 


0.075 


14 


60 


36 



As a consequence there is a significant cancellation of statistical errors when the 
force F(r) and the slope c(r) are calculated. It is then also profitable to keep track 
of the statistical correlations in the further analysis of the data. 



5. Simulation results 

In table 1 we list the values of the coupling f3 = 6a d ~ 4 /gQ and other important 
parameters of the simulations that we have performed. The entries in the last three 
columns are the maximal distance and the associated lattice sizes where we have 
accurate results for the Polyakov loop correlation function. At smaller distances, the 
lattices do not need to be as large, but were always chosen so that the systematic 
errors (other than the lattice effects) are sufficiently suppressed. Further algorithmic 
details and a collection of data tables are included in appendices A and B. 

In the following we first discuss the results at the smallest values of the lattice 
spacing and shall then argue (in subsect. 5.3) that the lattice effects at these points 
are already very small. 

5.1 F(r) and c(r) in four dimensions 

One of the striking results of our computations is that the force F(r) appears to 
be a linear function of 1/r 2 in the range r > 0.5 fm (see fig. 2). The errors on the 
data points shown in this figure increase from below 0.04% to about 0.12% at the 
left-most point, and the four points at the largest values of r lie on a straight line 
to this level of accuracy. 
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Fig. 2. Force F(r) in the four-dimensional theory versus 1/r 2 at (3 = 6.0. Physical 
units are set by the Sommer scale r = 0.5 fm [22] and the dotted line is a linear fit 
to the four points at the largest distances. Errors are smaller than the data symbols. 

There is still some curvature in the data, which is more clearly seen in fig. 3 where 
we plot the slope c(r). The errors on these points are at most 2.5% and significantly 
less than this at the smaller distances. Note that c(r) can be interpreted as a running 
coupling, 



c (r) = -| a (l/r), a(/x) 



m 2 

4vr ' 



(5.1) 



that satisfies the usual perturbative renormalization group equation 

r\_ OO 

^ n=0 



(5.2) 



The coefficients b n are known up to three-loop order in this scheme [29-33] , and using 
these we can integrate eq. (5.2) from r = 0, taking the A parameter determined in 
refs. [34,35] as the initial condition. 
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. 





0.5 1 r[fm] 

Fig. 3. Slope c(r) in three and four dimensions at f3 = 20.0 and (3 = 6.0 respectively. 
The curves on the left are obtained from perturbation theory, while those on the right 
derive from the string theory formula (3.9) with 6 = and b = 0.04 fm. 

In the related case of the force F(r), the range of applicability of the perturbation 
expansion extends to distances of around 0.2 or at most 0.3 fm [36]. We refer to this 
paper for further details and only remark that the error on the A parameter alone 
results in a significant spread of the three-loop curves shown in fig. 3 (shaded area 
on the left). It is in any case evident from the plot that c(r) breaks away from the 
perturbative behaviour near r = 0.2 fm. 

At larger distances the data show a shallow minimum and then gently increase to 
about —0.293 at r = 1 fm. This is still 12% away from 7 = —0.262, but not really 
inconsistent with the effective string theory, because the correction proportional 
to b in eq. (3.9) (and possibly higher-order terms) may not be negligible at these 
distances. As shown in fig. 3 we can in fact easily make contact with the data if we 
allow for such a correction. 
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5.2 Three-dimensional theory 

In three dimensions the static potential at short distances has a completely different 
shape, but the effective string theory predicts the same behaviour at large distances, 
with a coefficient 7 that is precisely half as big as in four dimensions. 

Our numerical results fully agree with this (see fig. 3). No b/r correction is required 
in this case to match the data points at the largest distances, while at the smaller 
values of r we observe a smooth transition from perturbative to string behaviour. 
In this theory the perturbation expansion to two-loop order reads [32] 

c(r) = -^g 2 r + 0(g 6 ), (5.3) 

where g denotes the gauge coupling in the continuum limit. To express the coupling 
in physical units, we used the conversion factor g 2 ro = 2.2 (which can be inferred 
from our data) and set ro = 0.5 fm as usual. 

5.3 Lattice effects 

In both three and four dimensions, we have simulated lattices at different values of 
the lattice spacing in order to study the approach of the force F(r) and the slope 
c(r) to the continuum limit (see table 1). The comparison of the data produced in 
these simulations requires all dimensionful quantities to be scaled to some physical 
units. It is natural to take the Sommer radius r as the basic reference scale in the 
present context, and we shall thus be interested in determining the lattice-spacing 
dependence of the dimensionless functions r 2 F(xr ) and c(xr ) in a range of the 
scale factor x. 

Since the distances xr are in general fractional multiples of the lattice spacing, 
an interpolation formula must be specified. Given the approximate shape of F(r) 
and c(r), we decided to use a three-point polynomial interpolation in 1/x 2 and 1/x 
respectively. Note that it would be logically incorrect to assign a systematic error to 
the interpolation, because the latter merely defines the functions between the lattice 
points. On the other hand, a poor choice of the interpolation formula will lead to 
artificially enhanced lattice effects. 

From the tables in appendix B it is clear that the lattice effects are quite small 
at the lattice spacings considered. In both three and four dimensions, the residual 
effects on the lattices with the smallest spacing are probably no more than about 
0.1% in the case of the force rgF(xro) and 1 — 2% in the case of the slope c(xr ). The 
theoretically expected scaling proportional to a 2 [37,4] is, however, not observed, at 
least not in four dimensions, where the lattice effects at (3 = 5.7 appear to be larger 
than what would be inferred from the other lattices. 
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For this reason, and since we have data at only three values of the lattice spacing, 
we do not attempt to perform an extrapolation to the continuum limit. We would 
in any case not expect the extrapolation to have any impact on the outcome of our 
work, because the lattice effects are small and because in four dimensions the values 
of c(xr ) at the larger distances xr would only be moved closer to 7. 



6. Conclusions 

The results reported in this paper show that string behaviour in the static potential 
sets in at quark-antiquark separations around 0.5 fm. In three dimensions the con- 
firmation of the string theory formula (1.1) is particularly impressive, while in four 
dimensions a small higher-order correction needs to be included to match the data 
points at the largest distances in fig. 3. 

The observed agreement with the effective string theory is somewhat surprising 
and in fact difficult to understand, because the physical picture of a thin fluctuating 
flux tube is hardly correct at these distances. Perhaps this is an indication of the 
existence of an exact dual formulation of SU(iV) gauge theories in terms of a funda- 
mental string theory, along the lines of ref. [38] for example, from which the effective 
theory derives. Note that not only the general form of the effective action, but also 
the values of the coefficients of the higher-dimensional terms would be predictable 
in this case (to the extent the fundamental string theory is tractable). 

Another logical possibility is that eq. (1.1) may have a different origin and that the 
proper interpretation of our results still needs to be found. It has been pointed out in 
this connection that the spectrum of the excited levels E n does not appear to agree 
with the one obtained at leading order of the effective string theory [25,26]. Whether 
this rules out the effective theory is not completely obvious, however, because the 
observed discrepancies (which are significant even at distances r > 2 fm) could be a 
result of the higher-order corrections. In particular, the accidental degeneracies of 
the spectrum at leading order probably disappear once the non-linear interactions 
are switched on. 

A systematic study of the higher-order effects would now be required to resolve 
this issue. Accurate calculations of the lowest energy values in each symmetry sector 
may then already provide enough information to determine the coefficients of the 
few leading terms in the effective action. It would evidently also be interesting to 
see whether our results carry over to other non-abelian gauge groups and to all 
representations of the Wilson lines with non-trivial transformation behaviour under 
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the centre of the gauge group. We finally note that the computation of c(r) at 
distances larger than reported here rapidly requires enormous amounts of computer 
time, because the significance loss in eq. (4.3) (which grows proportionally to or /a 2 ) 
must be compensated by increasing the statistics. 
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Appendix A. Simulation algorithm 

The multilevel algorithm that was used in the present study is described in detail 
in ref. [6]. Here we merely wish to add a few comments on performance and opti- 
mization issues. 

All tests of the algorithm reported in ref. [6] were carried out on four-dimensional 
lattices at f3 = 5.7. We now confirm that it performs well also at larger values of (3 
and in three dimensions. For illustration we quote the computation of the Polyakov 
loop correlation function at (3 = 6.0, T = 48a and r = 12a. This is an extreme 
case where the area of the minimal surface bounded by the loops reaches 5 fm 2 . The 
correlation function is very small at this point (about 1.1 x 10~ 25 ), but with one 
month of running time on a PC-cluster with 24 Pentium 4 processors it was possible 
to determine its value to an accuracy better than 3%. 

Contrary to what may be understood from ref. [6], the more complicated forms of 
the algorithm, with many levels of nested averages, are often not worth while. Most 
of our results have in fact been obtained with only one level, but a relatively large 
number (up to a few thousands) of "measurements" of the two-link operators on the 
associated time slices. An important side effect of such high numbers of time-slice 
updates is that the statistical correlations between the data at the distances covered 
in a given run tend to be enhanced (cf. subsect. 4.3). 

For future applications of the multilevel algorithm at larger distances r, it may be 
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necessary to study and improve its efficiency at the level of the time-slice averages. 
There is currently no understanding of the autocorrelation times in this subsystem, 
and rather than increasing the statistics, there may be better ways to cope with the 
exponential decay of the averaged two-link operators. 



Appendix B. Data tables 

In tables 2-7 we only list the most important results for the force F(r) and the 
slope c(r). More detailed data tables, including error correlation matrices, may be 
obtained from the authors in electronic form. 



Table 2. Simulation results at f3 = 6.0 in 4 dimensions 



rja 


rja 


a 2 F{f) 


rja 


c(f) 


3 


2.277 


0.102491(20) 


2.700 


0.28184(10) 


4 


3.312 


0.073854(22) 


3.729 


0.29724(18) 


5 


4.359 


0.062399(24) 


4.786 


0.30339(38) 


6 


5.393 


0.056871(25) 


5.833 


0.30515(70) 


7 


6.414 


0.053800(27) 


6.864 


0.3041(12) 


8 


7.428 


0.051924(30) 


7.886 


0.3033(17) 


9 


8.438 


0.050687(34) 


8.901 


0.3008(28) 


10 


9.445 


0.049834(38) 


9.912 


0.2992(80) 


11 


10.451 


0.049213(46) 






12 


11.455 


0.048839(72) 
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Table 3. Simulation results at ft — 20.0 in 3 dimensions 



r/a 


r/a 


a 2 F{f) 


r/a 


c(f) 


3 


2.379 


0.0508627(32) 


O OAO 

2.808 


0.081536(41 J 


4 


3.407 


0.0434942(32) 


3.838 


0.099772(77) 


5 


4.432 


0.0399646(33) 


4.875 


0.112270(96) 


6 


5.448 


0.0380274(35) 


5.902 


0.12057(13) 


7 


6.458 


0.0368546(38) 


6.920 


0.12596(23) 


8 


7.464 


0.0360941(41) 


7.932 


0.12945(41) 


9 


8.469 


0.0355749(45) 


8.941 


0.13100(55) 


10 


9.473 


0.0352080(49) 


9.948 


0.13184(69) 


11 


10.475 


0.0349401(55) 


10.953 


0.13276(92) 


12 


11.478 


0.0347382(62) 


11.957 


0.1321(17) 


13 


12.480 


0.0345838(73) 


12.961 


0.1271(34) 


14 


13.481 


0.0344673(92) 







Table 4. Values of rQF(xro) in 4 dimensions 



X 


= 5.70 


(5 = 5.85 


P = 6.00 


0.8 


1.82181(23) 


1.82203(57) 


1.82041(16) 


1.0 


1.65 


1.65 


1.65 


1.2 


1.55341(12) 


1.55602(34) 


1.55691(18) 


1.4 


1.49372(23) 


1.49983(49) 


1.50100(36) 


1.6 


1.45561(31) 


1.46369(60) 


1.46476(51) 


1.8 


1.42963(39) 


1.43920(76) 


1.44010(70) 


2.0 


1.41122(47) 


1.4217(10) 


1.4230(11) 
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Table 5. Values of r§F(xr ) in 3 dimensions 



X 


= 11.0 


= 15.0 


= 20.0 


0.8 


1.720701(26) 


1.719410(57) 


1.718083(52) 


1.0 


1.65 


1.65 


1.65 


1.2 


1.608316(37) 


1.610380(53) 


1.610988(73) 


1.4 


1.582508(55) 


1.58583(10) 


1.58689(12) 


1.6 


1.565497(77) 


1.56965(14) 


1.57112(16) 


1.8 


1.553860(90) 


1.55852(19) 


1.56022(22) 


2.0 


1.54555(11) 


1.55048(24) 


1.55271(34) 


Table 6. Values of c(xr ) in 4 dimensions 


X 


= 5.70 


= 5.85 


= 6.00 


1.0 


-0.31427(43) 


-0.3067(10) 


-0.30493(58) 


1.2 


-0.32150(57) 


-0.3048(14) 


-0.3045(11) 


1.4 


-0.32113(72) 


-0.3036(13) 


-0.3039(15) 


1.6 


-0.3176(10) 


-0.3003(22) 


-0.3014(24) 


1.8 


-0.3150(19) 


-0.2990(50) 


-0.2995(64) 


Table 7. Values of c(xro) in 3 dimensions 


X 


= 11.0 


0= 15.0 


= 20.0 


1.0 


-0.13225(14) 


-0.12696(18) 


-0.12505(21) 


1.2 


-0.13842(16) 


-0.13157(32) 


-0.12972(40) 


1.4 


-0.14106(30) 


-0.13475(40) 


-0.13144(55) 


1.6 


-0.14153(33) 


-0.13552(67) 


-0.13273(84) 


1.8 


-0.14106(52) 


-0.13654(97) 


-0.1316(18) 
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